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We investigate theoretically the ground-state configurations of two-dimensional charged-particle 
systems with an elliptical hard-wall boundary and their vibrational eigenmodes. The systems exhibit 
a series of structural transitions, finally changing from a zigzag structure to a one-dimensional 
Coulomb chain, as the eccentricity of the elliptical hard- wall boundary is increased. 
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I. INTRODUCTION 

In 1934, Wigner predicted that electrons will crystal- 
lize in a three-dimensional electron gas and form a lat- 
tice when the density of the gas is lowered to a certain 
critical value. Since then, Wigner crystals have been ob- 
served in various low-dimensional systems. These sys- 
tems include electrons on the surface of liquid helium |l| , 
electrons confined in semiconductor quantum dots 
colloidal suspensions with macroscopic particles vor- 
tices in superfluids 3 , and strongly coupled dusty plasma 
i 

Recently, classical Wigner crystals have aroused con- 
siderable interest both experimentally |E Q, H, H E| and 
theoretically [H El El Q El . Theoretical results for 
the ground-state configurations and normal mode spec- 
tra of classical Wigner crystals agree with experiments 
very well. Very recently, Melzer observed experimen- 
tally that a two-dimensional (2D) charged-particle 
system under an anisotropic parabolic confining poten- 
tial forms a one-dimensional (ID) chain through a zigzag 
transition, which agrees with the theoretical expectation 
[Til . Most previous studies have focused on systems 
under a parabolic confining potential or Coulomb poten- 
tial induced by a fixed positive charge. There have been 
only a few investigations of isotropic hard-wall confine- 
ment 0, 01 . These have adopted Monte Carlo simula- 
tion augmented with the Newton optimization technique. 
Since Newton optimization cannot directly minimize the 
energy of a system under the hard-wall confining po- 
tential, the authors E3l adopted an r"-like confinement 
to simulate hard-wall confinement in order to study the 
static and dynamical properties, because an r"-like po- 
tential approaches the hard- wall confining potential when 
the exponent n approaches the infinity. 

In this paper, we adop t a nonlinear optimization al- 
gorithm (L-BFGS-B) E3 to investigate the ground-state 
configuration of a 2D charged-particle system under an 
ellipitical hard-wall confining potential. The excitation 
modes for different shapes of confining potential are also 
studied. The present paper is organized as follows. In 
Sec. II, our model and numerical method are introduced. 



In Sec. Ill, we calculate the ground-state configurations 
and find that they exhibit a series of structural phase 
transitions as the anisotropy of the elliptical hard-wall 
boundary is increased. Sec. IV describes the influence 
of the anisotropy of the hard-wall confining potential on 
the normal mode excitation spectrum. Finally, our con- 
clusions are given in Sec. V. 



II. THEORETICAL MODEL 

The potential energy of a 2D system of N charged 
particles including the Coulomb potential under a hard- 
wall confining potential Vc{r) is given by 
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In the case of an elliptical hard-wall confining poten- 
tial, utilizing the coordinate transform Xi — Vi cos 9i, yi — 
evi sin 9i , where e refers to the eccentricity of the ellipiti- 
cal boundary, the elliptical hard-wall confining potential 
can be expressed as 



Vc{r) = { 



if r < i?c 
oo if r > i?c 
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where Rc is the major axis length. 

By choosing tq = Rc as the unit of length and Eq = 
/[sRc) as the unit of energy, the potential energy in 

dimensionless form is 
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where < r,, < 1, < OiOj < 2tt. From this form, 
we note that the Hamiltonian depends only on the num- 
ber of particles N and the eccentricity e of the ellipitical 
boundary. Due to the invariance of the system configu- 
ration after a 7r/2 rotation, we only investigate the cases 
of < e < 1 in this paper. 

The L-BFGS-B algorithm E3|) which is an im- 
proved version of the BFGS [1^ algorithm, can find 
the minimum of a function f[r{xi,X2, ■ ■ - Xi- • ■)] under 
the constraint of k < Xi < Ui. To solve the boundary 
problem, we obtain the piece- wise linear path x{t) = 



P (xi — tbi, I, u) by projecting the steepest descent direc- 
tion onto the feasible region: 

{k if Xi < k 
Xi if Xi e [Zi, Ui] . (4) 
Ui if Xi > Ui 

Such a transformation enables us to find the local 
minimum of the univariate piece- wise quadratic Qnit) = 
yn(x{t)), in which x'^ is defined as the Cauchy point. 
The variables whose value at x'^ is at the lower or upper 
bound, comprising the active set A{x'^)^ are held fixed. 
Then we consider the following problem: 

min{f{x : Xi = x^, Vi £ A{x'^)} (5) 
k<x,< u, yi ^ A(a;=). (6) 

First, adopting the conjugate gradient algorithm 
on the subspace of free variables, we approximately solve 
Eq. ISJl. In the iterative process, we take as the initial 
value for iteration, and then truncate the path to satisfy 
Eq. ((BJ. After an approximate solution fn+i has been 
obtained, / decreases further along fi„ — r„+i — r„ on the 
condition that r„_|_i is in the appropriate range (bound- 
ary constraint) . Then the approximation Hessian matrix 
and the gradient vector at r„+i are calculated, and a 
new iteration starts. Such iterations continue until the 
required accuracy is achieved. More detailed informaton 
about this algorithm can be found in Ref. jl7|. 

The total energy of the system under the elliptical 
hard-wall confining potential (see Eq. JSJ ) can be min- 
imized by the L-BFGS-B algorithm. First, all particle 
coordinates are generated randomly. Secondly, the L- 
BFGS-B algorithm is adopted to minimize the total en- 
ergy of the system. We repeat the above procedure many 
times. Finally, we choose the configuration with the low- 
est energy as the ground-state of the system. To check 
the accuracy and validity of this algorithm, we repro- 
duced the ground-state configurations for an isotropic 
hard- wall potential TT] . 

Once the system ground-state is obtained, the vibra- 
tional eigenmode can be calculated by diagonalizing the 
dynamical matrix (see Appendix). The eigenfrequencies 
in this paper are expressed in the unit ujq = \/ Eo/{mrQ). 

III. NUMERICAL RESULTS AND 
DISCUSSIONS 

A. The ground-state configuration 

When the charged particles are confined by a hard-wall 
potential, they repel each other due to the Coulomb in- 
teraction so that some particles first move to the bound- 
ary until the outer shell is fully occupied, the remaining 
particles stay in the inner shells. The shape of the bound- 
ary significantly affects the system ground-state config- 
uration. Fig.l shows the evolution of the ground-state 
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FIG. 1: The ground state configurations of three clusters with 
difi'erent numbers of particles — 12, 26, 50 fors different 
eccentricities of the elliptical boundary. 



configurations for different numbers of particles as the 
boundary eccentricity e decreases. The ground-state con- 
figurations in the anisotropic case are quite different from 
those in the isotropic case (e = 1). They change signifi- 
cantly with decreasing eccentricity e, and exhibit zigzag 
structures when e approaches zero. 

Taking the iV = 50 cluster as an example, we can find a 
continuous deformation from circular to elliptical rings as 
the eccentricity e decreases. The inner ring collapses into 
a line with decreasing e. Notice that more and more par- 
ticles leave the central line and merge into the adjacent 
outer ring. Then the system forms two elliptical rings 
(see Fig. 1). This process repeats until the system ex- 
hibits a quasi-one-dimensional zigzag configuration and 
becomes a straight line at e = 0. In addition, inner shells 
are more sensitive to the change of the eccentricity e than 
outer shells. Decreasing e further, particles in the inner 
shell move to the outer shell, and the inner shell finally 
collapses into a line. 

We find that there are two different kinds of transi- 
tions: slow and abrupt transitions. For convenience, we 
refer to the configuration in which — 1 particles are on 
the outer shell and one at the centre as a-type, and the 
configuration in which all particles are on the outer shell 
as /3-type. From the Table I in Ref. we know that 
the system exhibits the a- type configuration for A^ = 12 
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FIG. 2: The ground-state energys , the first derivative and 
second derivative of the ground state energy with respect to 
the eccentricity as a function of the eccentricity for A'^ = 12. 
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FIG. 4: Same as Fig. 3, but for A*" — 50. The numbers denote 
the number of particles in the inner shell. 
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FIG. 3: The second-order derivative d^E'/de'^ as a function 
of eccentricity for A'^ = 26, where b-h correspond to the con- 
figurations for A*' = 26 in Fig.l 
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FIG. 5: Density per unit length n{x), and zigzag angle as a 
funtion of coordinate x for different eccentricities. 



under a circular hard-wall confining potential. But for an 
ellipitical confining potential, the transition from a-type 
to /3-type occurs as the boundary eccentricity e decreases 
to a certain value. Fig. 2(a) shows the system energy 
of TV = 12 as a function of the eccentricity e when the 
major axis is fixed and the minor axis decreases with e. 
When e decreases and the area is compressed, the total 
energy increases. We adopt E' = E ■ e^/^ and calculate 
the first and second derivatives of E' with respect to e, 
i.e. (3ii^'/9e, 9^ i?'/9^e, which are shown in Fig. 2(b). An 
obvious dip in d^E'/d^e appears with the decrease of e. 
As seen in Fig. 1, the system ground-state is a-type at 
e = 0.7 and /3-type at e = 0.65. The transition from a- 
type to /?-type takes place at e w 0.666, 0.546, and 0.483 
respectively for N = 12, 13, and 14. 

In Figs. 3 and 4, we plot the second derivative of en- 



ergy with respect to the eccentricity e for iV = 26 and 50 
as a function of eccentricity e. More transitions, which 
are denoted by the vertical lines, can be found. From the 
evolution of three clusters, the ground-state configuration 
finally transits from a-type to /?-type and finally to a ID 
chain through a zigzag structure. Note how the system 
evolves from 2D to ID. In Fig. 5(a) and 5(b), we plot 
zigzag angles and particle density of zigzag structures and 
a ID chain as a function of particle position for A^ = 50. 
From Fig. 5(a), we find the density per unit length of 
a ID chain (e = 0) exhibits behavior similar to that of 
the zigzag structure as the x coordinate increases. It re- 
mains almost constant in the central region and increases 
significantly at the edges of the chain. These behaviors 
of a ID chain under hard- wall confinement are very dif- 
ferent from those under a parabolic potential, which was 
investigated by Dubin. |^ and Melzer Q. 
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FIG. 6: The frequency spectra as functions of the eccentric- 
ity e for a cluster N = 3. The dashed, soUd, and dotted Unes 
denote the frequencies of three vibrational modes. The eigen- 
vectors, with frequencies arranged in ascending order from 
left to right, arc in the insets in Fig. 6(a) for e — 0.1. The 
eigenvectors of e = 1.6 are in the lower insets of Fig. 6(b), 
and those of e = 2.8 are in the upper insets of Fig. 6(b). 



IV. THE EIGENMODE SPECTRUM 

In this section, we discuss the vibrational modes of 
the ground state under the elliptical hard-wall confining 
potential. 

Because the radial motion of the particles at the 
boundary is practically frozen, only 2N — Ne modes exist 
when the Ne particles are located at the boundary. But 
under ellipse hard-wall confinement, the rotational mode 
with zero frequency at e = becomes nonzero and some 
degenerate modes split due to the breaking of cylindri- 
cal symmetry. In addition, increasing the eccentricity e 
results in abrupt transitions of the ground-state config- 
urations {N > 12), and their eigenmodes change signif- 
icantly because of the abrupt change of the number of 
particles at the boundary. 
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FIG. 7: The frequency spectra as funtions of the eccentricity 
e for a cluster ISS = 12. 



nism becomes more important, and two split HF modes 
degenerate approximately again and the low-frequency 
(LF) mode saturates gradually. From the eigenvectors 
at e = 0.1 (see the insets of Fig. 6(a)), we find that 
decreased eccentricity e enhances the confinement along 
the y axis, and consequently leads to the increases of fre- 
quencies of the vibration modes along the y axis. Thus 
those frequencies tend to infinity in the extreme case of 
e = 0. From the ID Coulomb chain potential fimction, 
we diagonalize the dynamical matrix of the ground-state 
configuration for an A'' = 3 cluster and get that the lowest 
frequency is 2.0 and the two high frequencies are infinity, 
which agrees with Fig. 6(a). In addition. Figs. 6(a) and 
6(b) show that the two split modes cross respectively at 
e « 2.3 and 0.44 (notice 2.3 « 1/0.44). We compare the 
eigenvector of the highest frequency at e = 1.6 with that 
at e = 2.8 (see the insets of Fig. 6(b)). One can see that 
two particles near the major axis move in the same direc- 
tion at e = 1.6 but in the opposite direction at e = 2.8. 
This verifies that the cross occurs indeed between e = 1.6 
and e = 2.8. 



A. N = 3 pEirticIe cluster 

The frequencies of three vibrational modes in an A'' = 3 

cluster as a function of eccentricity e arc shown in Figs. 
6(a) and 6(b). One can find two different mechanisms 
driving the transition of the ground state configuration 
as the anisotropy increases. One is due to anisotropic 
hard-wall confinement. The other is the minor (major) 
axis of the ellipse decreases (increases), and consequently 
the Coulomb interaction increases (decreases) for e < 1 
(e > 1). When e approaches 1, the former mechanism 
is dominant and the rotational mode frequency increases 
(see dashed line in Fig. 6(b)). In addition, this mech- 
anism make high-frequency (HF) degenerate mode split 
into two modes (see Fig.6(a,b) solid and dotted lines). 
When e decreases (increases) further, the latter mecha- 



B. N = 12 particle cluster 

In this subsection we discuss a larger system {N = 12) 

to illuminate the effect of the elliptical hard-wall bound- 
ary on the spectrum of the system. Fig. 7 shows all 
eigenfrequencies of the iV = 12 cluster as a function of 
the eccentricity e. From this figure, we see that the fre- 
quency spectra exhibit many jumps between e = 0.668 
to e = 0.6. This is because the system transits from a- 
type to /3-type at e « 0.666 (see the definition of a-type 
and /3-type in Section II). The distance between neigh- 
boring particles decreases, resulting in an abrupt increase 
of the eigenfrequency. In addition, the highest frequency 
mode diminishes due to the freezing of the radial mo- 
tion of all particles at the boundary (see the thirteenth 
mode in Fig. 7(b)). It is interesting to note that a small 
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dip appears at e « 0.22 which cannot be found in the 
energy figure of the system (see Fig. 2). Checking the 
ground-state configuration around e = 0.22, we find that 
the system changes from /3-type to the zigzag structure 
(see the insets in Fig. 7(a)). Next, we explain the effect 
of the eUiptical boundary on vibrational modes in detail 
in three regions. 

When e > 0.66, the system ground state configuration 
is a-type, and there are thirteen vibration modes. Among 
these modes, the anisotropic effect results in the apparent 
splitting of the LF modes, but the effect on the other 
higher frequency modes is negligible. 

When 0.22 ^ e < 0.66, the system ground state con- 
figuration is /3-type. Because the number of particles 
increases by one at the boundary, and the HF modes re- 
compose due to the configuration transition from a-type 
to /3-type. As the eccentricity e decreases, the second 
mode crosses with the third mode (see Fig. 7(a)). A simi- 
lar feature was observed experimentally in an anisotropic 
parabohc confinement 0. In the HF part (see Fig. 7(b)), 
the ninth mode departs from the eighth mode and be- 
comes degenerate with the tenth mode gradually, corre- 
sponding to the slow transition of the configuration. 

When e < 0.22, the system ground-state configuration 
is a zigzag structure. From the inset in Fig. 7, we find 
that the zigzag structure destroys the configuration sym- 
metry, resulting in the splitting of the degenerate modes 
in the isotropic hard-wall confining potential. The low- 
est vibrational mode is no longer the rotation mode, and 
more crosses between LF modes (see Fig. 7(a)) can be 
found. The LF eigenvectors at e = 0.1 and e = 0.5 show 
clearly the crosses of different modes. From the eigen- 
vectors of the = 12 cluster in Fig. 8, we can compare 
the vibration direction of eigenmodes at e — 0.1 with 
that at e = 0.5 and find that all corresponding particles 
move in the same way among the vibration modes (a) 
and (o), (b)and (q), (c) and (m), (d) and (r). The par- 
ticles with larger vibration amplitudes in the LF modes 



appear close to the center of the system as e decreases. 
But the behavior of the HF modes is opposite to that of 
the LF modes (see Fig. 8). 



V. CONCLUSION 

The ground-state configuration and dynamical prop- 
erties of a 2D classical cluster of charged particles in an 
elliptical hard-wall confining potential are investigated 
theoretically. The dependence of the ground-state con- 
figuration on the number of particles and the eccentricity 
of the confining potential are discussed in detail. As the 
eccentricity e decreases, the system exhibits two differ- 
ent kinds of transitions, i.e., slow and abrupt transitions. 
When increasing the system anisotropy, circular shells 
are deformed into ellipses continuously, the inner shell 
collapses into a line, and then number of particles in 
the line decreases gradually. Last, the system exhibits 
a zigzag structure and forms a ID Coulomb chain when 
e approaches zero. 

The boundary anisotropy significantly influence the 
frequency of the vibration modes. With increasing the 
anisotropy of the boundary, the rotation mode becomes 
nonzero and degenerate modes split due to the break- 
ing of cylindrical symmetry. The rotation mode is no 
longer the lowest vibration mode when the eccentricity e 
approaches zero. 
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In this appendix, we present the dynamical matrix 
of a 2D classical system under the elliptical coordi- 
nates. Using the coordinate transformation, Xi — 
epi cos di, Hi — Pi sin 0i, where e denotes the eccentric- 
ity of the elliptical hard-wall boundary, one can obtain 
X — epcosO — {psm6)9 and y = psm9 + {p cos 9)9. 
The equilibrium positions of the particles are denoted by 

In the harmonic approximation around the equi- 
librium position, the Lagrangian can be expressed 



as L({p,},{0J) = 



e cos 



sm 



pf{e'^sm^9i + cos'^9t 
^Ei'^P.,P,Sp^Sp, 



2pi{l - e^) sin 9. 1 cos 9 ip9) 



2$„ 



5 Pi = 



Pi 



- Pi 



5pi59j 



Pi.Pj 



Since S9i, 5pi 



59iS9j), where 

^dp^dpj ''e?' 

0, linearizing 



(-21^\ Pic 

the equation of motion for {pi} and {9i} and defining 
T = {5p, p^59), we can obtain the equation of the par- 

,2 

tides' movement under matrix form = —M ■ T. 

So the dynamical matrix of a 2D classical system Af 

B G C-E B H C-F _ 

where 



IS 



AC-BD 

A-G 
BD-AC 



AC~BD 
BE 



AC-BD 
A-H 



B-D-A-C BD-AC 



AC-BD 

BF 
BD-AC 
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K,] = + [B.,.] = {1 - e') sm 6^ cos 0^ , [/f. .] = 

[Ci,i] = (e2 sin2 6*^ +cos2 0i), = (1 - e^) sin cos 61°, 

[E^J] = -$p.,p„ [^^ij] = [G,,] = and 
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